Getting Started


Setting up the python environment

python has been chosen as the main language of this book, for several reasons. It doesn’t have too much jargon, and the signal-to-noise ratio in code is very high (it’s almost pseudocode). While python is not the optimal language for scientific computing, it is familiar to most people with basic knowledge of programming, and the learning curve for python is not as steep as C++, Java, or Julia. One of the drawbacks associated with using python is that it doesn’t handle large graphs or datasets very well. C++ and Java both deal with larger structures more efficienly, and many authors of seminal papers in the field of GIS provide C++ and Java implementations alongside their papers. Julia, on the other hand, is a perfect midpoint between python and C++/Java. It maintains the readability of python but was developed specifically for scientific computing applications.

For the sake of accessibility, the examples and exercises in this book will be based entirely in python. Implementations in other languages are welcome as pull requests. Simply make a PR of this repo and we’ll include any submissions that look promising.

From this point forward, the book will assume that you already have python 3.6+ installed on your system. For installation instructions specific to your operating system, see this Beginner’s Guide.

There are two primary ways to install the required packages you’ll need: pip and conda. It is recommended to use pip as conda may have some adverse effects on system dependencies when used improperly in Linux.

Note

If you intend on developing in Windows, it is also recommended that you leverage the convenience of the Windows Subsystem for Linux (WSL). This allows you to use the full capabilities of Linux from within Windows. You can read more about setting up WSL here. After enabling WSL, you can proceed with the rest of this book as if you were operating in Linux.

Using pip3

In the terminal, execute the following commands:

$ sudo apt update
$ sudo apt install python3-pip

Install venv and create a python virtual environment:

$ sudo apt install python3.8-venv
$ mkdir <new directory for venv>
$ python3 -m venv <path to venv directory>

Make sure that you replace python3.8 with the version of python you are using.

You can now access your venv using the following command:

$ source <path to venv>/bin/activate

Using the MacOS terminal:

$ python3 -m ensurepip --upgrade

venv is included with python 3.8+. You can run the following command to create a virtual environment:

$ mkdir <new directory>
$ python3 -m venv <path to venv directory>

You can now access your venv using the following command:

$ source <path to venv>/bin/activate

pip should come preinstalled with Windows installations of python for versions newer than 3.4.

Make sure to check the “Add Python to PATH” option when using the Python installer. If you don’t, you’ll need to add python, pip, and various other programs to the system path manually.

After installing python, you may need to restart your computer in order for the changes to take effect.

To create a virtual environment, execute the following either in command prompt, or Windows Powershell:

C:> py -m venv <new directory>

To activate the venv:

C:> cd <your venv directory>
C:> .\Scripts\activate

Note

There are several things to note when using python in Windows.

  1. Using the py or python command may open the Windows Store app in newer releases of Windows. This can be disabled by going to “Apps and Features”, selecting “Application Execution Aliases”, and disabling the sliders for any application related to python (i.e. python3.8.exe, py.exe).

  2. pip3 is not sometimes not added to the system PATH by default. You may choose to add it manually, or simply use pip, as it should link to the same program.

  3. You can use python by calling the py command. If you prefer using python or python3, you may need to add these to the system PATH manually.

Using conda

Install conda for your OS using the guide found here.

You can create a conda environment like this:

$ conda create --name <name of env> python=<your version of python>

Access the newly-created environment using this command:

$ conda activate <your env name>

Installing Jupyter Notebook

All of the code in this book is stored in Jupyter Notebooks (.ipynb). To access these files directly, you have two options:

  1. Open Binder or Google Colab using the icon on the top right of a page containing python code. You can try it out with the current page.

  2. Install Jupyter Notebook locally.

Jupyter Notebook can be installed as follows:

$ pip3 install jupyterlab
$ pip3 install notebook
$ conda install -c conda-forge jupyterlab
$ conda install -c conda-forge notebook

Getting the data

Most of the data used in this book will be sourced from OpenStreetMaps. For any other datasets that are not OSM-related, you can download the data in whatever format it is available in, and import it into python using the appropriate methods. Open Data websites hosted by various government bodies are a great source of data related to infrastructure, population metrics, and land use. See the Datasets page for some examples.

For OpenStreetMaps, there are two primary ways of retrieving the data:

  1. Download the data as-is from OpenStreetMaps’ website and use tools like osmfilter to tune it as needed. This is not recommended as it is more difficult and not very efficient.

  2. Use OpenStreetMaps’ API (Overpass API) to query for data. This filters the data on the fly and you only retrieve what you need. The API is accessible in both Java and python.

Data Completeness

The data from OSM is not always “complete”. This doesn’t mean that there are major uncharted regions, but rather that neighbouring nodes are not always grouped correctly. For some nodes where there are feasible connections between them in real life, OSM still represents them as having no way or relation connecting them. This means that using the osmnx parser will result in the nodes being placed in separate graph components, which is not accurate to real-world conditions. We can use osrm to find routes between these kind of nodes and thus “complete” the graphs.

You can read more about the completeness of OpenStreetMaps data here:

  1. Completeness

  2. Completeness Metrics

The data model of OpenStreetMaps is surprisingly simple and consists only of three elements:

  1. Node represents a specific point on the earth’s surface defined by its latitude and longitude

  2. Way is an ordered list of between 2 and 2,000 nodes that define a polyline. Ways are used to represent linear features such as rivers and roads

  3. Relation is a multi-purpose data structure that documents a relationship between two or more data elements (nodes, ways, and/or other relations). Examples include:

    • A route relation, which lists the ways that form a major (numbered) highway, a cycle route, or a bus route.

    • A turn restriction that describes if a turn can be made from one way onto another.

    • A multipolygon that describes an area (whose boundary is the ‘outer way’) with holes (the ‘inner ways’).

All of the above can be found easily on the linked elements page, but there are two things you should be aware of:

  1. All ID’s of the same element type are unique globally, but they are not unique across element types (you can find a Way with the same ID as a Node).

  2. Ways and Relations are made by listing and referring to the ID’s of the Nodes that constitute them.


Basic Operations

This section contains a non-exhaustive list of operations on geospatial data that you should familiarize yourself with. More information can be found by consulting the Tools and Python Libraries page or the respective libary’s API documentation.

Creating a Graph from a named place

osmnx can convert a text descriptor of a place into a networkx graph. Let’s use the University of Toronto as an example:

import osmnx

place_name = "University of Toronto"


# networkx graph of the named place
graph = osmnx.graph_from_address(place_name)

# Plot the graphs
osmnx.plot_graph(graph,figsize=(10,10))
../_images/GettingStarted_7_0.png
(<Figure size 720x720 with 1 Axes>, <AxesSubplot:>)

The graph shows edges and nodes of the road network surrouding the University of Toronto’s St. George (downtown Toronto) campus. While it may look visually interesting, it extends a bit too far off campus, and lacks the context of the street names and other geographic features. Let’s restrict the scope of the network to 500 meters around the university, and use a folium map as a baselayer.

graph = osmnx.graph_from_address(place_name, dist=500)
osmnx.folium.plot_graph_folium(graph)
Make this Notebook Trusted to load map: File -> Trust Notebook

Suppose you want to get from one location to another on this campus.

Starting at the King Edward VII Equestrian Statue near Queen’s Park, you need to cut across the campus to attend a lecture at the Bahen Centre for Information Technology. Later in this section, you will see how you can calculate the shortest path between these two points.

For now, let’s plot these two locations on the map using ipyleaflet:

from ipyleaflet import *

# UofT main building
center=(43.662643, -79.395689) 
# King Edward VII Equestrian Statue
source_point = (43.664527, -79.392442)  
# Bahen Centre for Information Technology at UofT
destination_point = (43.659659, -79.397669) 

m = Map(center=center, zoom=15)
m.add_layer(Marker(location=source_point, icon=AwesomeIcon(name='camera', marker_color='red')))
m.add_layer(Marker(location=center, icon=AwesomeIcon(name='graduation-cap')))
m.add_layer(Marker(location=destination_point, icon=AwesomeIcon(name='university',marker_color='green')))
m

Notice that the way we create maps in ipyleaflet is different from folium. For the latter, the code is as follows:

import folium
m = folium.Map(location=center, zoom_start=15)
folium.Marker(location=source_point,icon=folium.Icon(color='red',icon='camera', prefix='fa')).add_to(m)
folium.Marker(location=center,icon=folium.Icon(color='blue',icon='graduation-cap', prefix='fa')).add_to(m)
folium.Marker(location=destination_point,icon=folium.Icon(color='green',icon='university', prefix='fa')).add_to(m)
m
Make this Notebook Trusted to load map: File -> Trust Notebook

Extracting Graph Information

Various properties of the graph can be examined, such as the graph type, edge (road) types, CRS projection, etc.

Graph Type

type(graph)
networkx.classes.multidigraph.MultiDiGraph

Note

You can read more about MultiDiGraph here.

Edges and Nodes

We can extract the nodes and edges of the graph as separate structures.

nodes, edges = osmnx.graph_to_gdfs(graph)

nodes.head(5)
y x highway street_count geometry
osmid
20964579 43.667513 -79.399806 traffic_signals 4 POINT (-79.39981 43.66751)
20979738 43.665937 -79.391876 traffic_signals 4 POINT (-79.39188 43.66594)
20979740 43.666067 -79.392493 NaN 3 POINT (-79.39249 43.66607)
20979742 43.666368 -79.393159 NaN 3 POINT (-79.39316 43.66637)
20979743 43.665877 -79.393329 NaN 3 POINT (-79.39333 43.66588)
edges.head(5)
osmid lanes name highway maxspeed oneway length geometry access service tunnel bridge width
u v key
20964579 389677900 0 3998177 3 Bloor Street West primary 50 False 8.362 LINESTRING (-79.39981 43.66751, -79.39971 43.6... NaN NaN NaN NaN NaN
390548868 0 4212261 2 St George Street tertiary 30 False 7.306 LINESTRING (-79.39981 43.66751, -79.39978 43.6... NaN NaN NaN NaN NaN
389678082 0 5090496 2 St George Street residential 40 False 12.405 LINESTRING (-79.39981 43.66751, -79.39986 43.6... NaN NaN NaN NaN NaN
389677899 0 234365457 2 Bloor Street West primary 50 False 9.046 LINESTRING (-79.39981 43.66751, -79.39991 43.6... NaN NaN NaN NaN NaN
20979738 2144928236 0 7685256 NaN NaN footway NaN False 9.967 LINESTRING (-79.39188 43.66594, -79.39182 43.6... NaN NaN NaN NaN NaN

We can further drill down to examine each individual node or edge.

# Rendering the 2nd node
list(graph.nodes(data=True))[1]
(20979738,
 {'y': 43.6659368,
  'x': -79.3918763,
  'highway': 'traffic_signals',
  'street_count': 4})
# Rendering the 1st edge
list(graph.edges(data=True))[0]
(20964579,
 389677900,
 {'osmid': 3998177,
  'lanes': '3',
  'name': 'Bloor Street West',
  'highway': 'primary',
  'maxspeed': '50',
  'oneway': False,
  'length': 8.362})

Street Types

Street types can also be retrieved for the graph:

print(edges['highway'].value_counts())
footway                       2030
service                        418
residential                    186
tertiary                       107
path                            82
secondary                       79
[steps, footway]                58
pedestrian                      36
unclassified                    24
primary                         22
steps                           20
[footway, service]              12
secondary_link                  10
cycleway                         6
[footway, path]                  4
[steps, path]                    4
[steps, service, path]           4
[steps, footway, corridor]       2
Name: highway, dtype: int64

GeoDataFrame to MultiDiGraph

GeoDataFrames can be easily converted back to MultiDiGraphs by using osmnx.graph_from_gdfs.

new_graph = osmnx.graph_from_gdfs(nodes,edges)
osmnx.plot_graph(new_graph,figsize=(10,10))
../_images/GettingStarted_25_0.png
(<Figure size 720x720 with 1 Axes>, <AxesSubplot:>)

Calculating Network Statistics

We can see some basic statistics of our graph using osmnx.basic_stats.

osmnx.basic_stats(graph)
{'n': 1090,
 'm': 3104,
 'k_avg': 5.695412844036698,
 'edge_length_total': 92073.89100000016,
 'edge_length_avg': 29.6629803479382,
 'streets_per_node_avg': 3.053211009174312,
 'streets_per_node_counts': {0: 0, 1: 154, 2: 0, 3: 579, 4: 349, 5: 7, 6: 1},
 'streets_per_node_proportions': {0: 0.0,
  1: 0.14128440366972478,
  2: 0.0,
  3: 0.5311926605504587,
  4: 0.3201834862385321,
  5: 0.006422018348623854,
  6: 0.0009174311926605505},
 'intersection_count': 936,
 'street_length_total': 48561.08599999995,
 'street_segment_count': 1613,
 'street_length_avg': 30.106066955982612,
 'circuity_avg': 1.0297899024310946,
 'self_loop_proportion': 0.0012399256044637321}

We can also see the circuity average. Circuity average is the sum of edge lengths divided by the sum of straight line distances. It produces a metric > 1 that indicates how “direct” the network is (i.e. how much more distance is required when travelling via the graph as opposed to walking in a straight line).

# osmnx expects an undirected graph
undir = graph.to_undirected()
osmnx.stats.circuity_avg(undir)
1.0328466941401129

Extended and Density Stats

Density-based statistics requires knowing the area of the graph. To do this, the convex hull of the graph is required.

convex_hull = edges.unary_union.convex_hull
convex_hull
../_images/GettingStarted_31_0.svg
import pandas
area = convex_hull.area
stats = osmnx.basic_stats(graph, area=area)
extended_stats = osmnx.extended_stats(graph,ecc=True,cc=True)
stats.update(extended_stats)

# Show the data in a pandas Series for better formatting
pandas.Series(stats)
/home/yinan/book-env/lib/python3.8/site-packages/osmnx/stats.py:405: UserWarning: The extended_stats function has been deprecated and will be removed in a future release. Use NetworkX directly for extended topological measures.
  warnings.warn(msg)
n                                                                                   1090
m                                                                                   3104
k_avg                                                                           5.695413
edge_length_total                                                              92073.891
edge_length_avg                                                                 29.66298
streets_per_node_avg                                                            3.053211
streets_per_node_counts                 {0: 0, 1: 154, 2: 0, 3: 579, 4: 349, 5: 7, 6: 1}
streets_per_node_proportions           {0: 0.0, 1: 0.14128440366972478, 2: 0.0, 3: 0....
intersection_count                                                                   936
street_length_total                                                            48561.086
street_segment_count                                                                1613
street_length_avg                                                              30.106067
circuity_avg                                                                     1.02979
self_loop_proportion                                                             0.00124
node_density_km                                                    10289596259367.783203
intersection_density_km                                             8835836787860.775391
edge_density_km                                                       869177215063338.25
street_density_km                                                      458416485189391.5
avg_neighbor_degree                    {20964579: 3.75, 20979738: 3.0, 20979740: 2.5,...
avg_neighbor_degree_avg                                                         3.126391
avg_weighted_neighbor_degree           {20964579: 0.4041057140547967, 20979738: 0.129...
avg_weighted_neighbor_degree_avg                                                 0.18964
degree_centrality                      {20964579: 0.0073461891643709825, 20979738: 0....
degree_centrality_avg                                                            0.00523
clustering_coefficient                 {20964579: 0, 20979738: 0, 20979740: 0, 209797...
clustering_coefficient_avg                                                      0.041101
clustering_coefficient_weighted        {20964579: 0, 20979738: 0, 20979740: 0, 209797...
clustering_coefficient_weighted_avg                                             0.005141
pagerank                               {20964579: 0.0005297849444390717, 20979738: 0....
pagerank_max_node                                                              390545025
pagerank_max                                                                    0.003326
pagerank_min_node                                                              306721057
pagerank_min                                                                    0.000138
eccentricity                           {20964579: 1391.2359999999999, 20979738: 1330....
diameter                                                                        1744.009
radius                                                                           878.559
center                                                                      [2143404200]
periphery                                                                   [2143499057]
closeness_centrality                   {20964579: 0.0013961273083024036, 20979738: 0....
closeness_centrality_avg                                                         0.00153
dtype: object

Warning

extended_stats has already been deprecated from osmnx. In the absence of an alternative, you can access these metrics directly through networkx.

CRS Projection

You can also look at the projection of the graph. To find out more about projections, check out this section. Additionally, you can also reproject the graph to a different CRS.

edges.crs
<Geographic 2D CRS: EPSG:4326>
Name: WGS 84
Axis Info [ellipsoidal]:
- Lat[north]: Geodetic latitude (degree)
- Lon[east]: Geodetic longitude (degree)
Area of Use:
- name: World.
- bounds: (-180.0, -90.0, 180.0, 90.0)
Datum: World Geodetic System 1984 ensemble
- Ellipsoid: WGS 84
- Prime Meridian: Greenwich
merc_edges = edges.to_crs(epsg=3857)
merc_edges.crs
<Projected CRS: EPSG:3857>
Name: WGS 84 / Pseudo-Mercator
Axis Info [cartesian]:
- X[east]: Easting (metre)
- Y[north]: Northing (metre)
Area of Use:
- name: World between 85.06°S and 85.06°N.
- bounds: (-180.0, -85.06, 180.0, 85.06)
Coordinate Operation:
- name: Popular Visualisation Pseudo-Mercator
- method: Popular Visualisation Pseudo Mercator
Datum: World Geodetic System 1984 ensemble
- Ellipsoid: WGS 84
- Prime Meridian: Greenwich

Shortest Path Analysis

Let’s revisit our trip across campus from the statue in Queen’s Park to our lecture hall at the Bahen centre.

To calculate the shortest path, we first need to find the closest nodes on the network to our starting and ending locations.

import geopandas

X = [source_point[1], destination_point[1]]
Y = [source_point[0], destination_point[0]]
closest_nodes = osmnx.distance.nearest_nodes(graph,X,Y)

# Get the rows from the Node GeoDataFrame
closest_rows = nodes.loc[closest_nodes]

# Put the two nodes into a GeoDataFrame
od_nodes = geopandas.GeoDataFrame(closest_rows, geometry='geometry', crs=nodes.crs)
od_nodes
y x highway street_count geometry
osmid
1907446268 43.664507 -79.392304 NaN 3 POINT (-79.39230 43.66451)
1633421938 43.659750 -79.398047 NaN 1 POINT (-79.39805 43.65975)

Let’s find and plot the shortest route now!

import networkx

shortest_route = networkx.shortest_path(G=graph,source=closest_nodes[0],target=closest_nodes[1], weight='length')
print(shortest_route)
[1907446268, 55808205, 55808194, 1907446267, 8699033082, 8699033084, 6542457312, 4953810914, 55808233, 299625330, 389677953, 7967019556, 7967019566, 4923076695, 55808571, 55808582, 389678112, 389678113, 389678146, 2143434862, 2143434860, 7311083158, 1258707987, 389678121, 50885147, 389678122, 389677906, 50885141, 389678180, 2143436415, 2143494216, 2143494214, 389678185, 1633421950, 389678184, 389678183, 389678216, 389678215, 389678226, 1633421933, 1633421938]
osmnx.plot_graph_route(graph,shortest_route,figsize=(15,15))
../_images/GettingStarted_41_0.png
(<Figure size 1080x1080 with 1 Axes>, <AxesSubplot:>)

Now let’s try visualizing this on a leaflet map. We’ll use the ipyleaflet library for this. Graphs rendered in ipyleaflet tend to be very slow when there are many nodes. For this reason, the code will fall back to folium if there are more than 1,500 nodes in the graph (this particular graph should have ~1,090 nodes).

Let’s make a map that shows the above route, with both starting and ending nodes shows as markers.

import ipyleaflet
import networkx

def draw_route(G, route, zoom=15):
    if len(G)>=1500:
        # Too many nodes to use ipyleaflet
        return osmnx.plot_route_folium(G=G,route=route)
    
    # Calculate the center node of the graph
    # This used to be accessible from extended_stats, but now we have to do it manually with a undirected graph in networkx
    undir = G.to_undirected()
    length_func = networkx.single_source_dijkstra_path_length
    sp = {source: dict(length_func(undir, source, weight="length")) for source in G.nodes}
    eccentricity = networkx.eccentricity(undir,sp=sp)
    center = networkx.center(undir,e=eccentricity)[0]

    # Create the leaflet map and center it around the center of our graph
    G_gdfs = osmnx.graph_to_gdfs(G)
    nodes_frame = G_gdfs[0]
    ways_frame = G_gdfs[1]
    center_node = nodes_frame.loc[center]
    location = [center_node.y,center_node.x]
    m = ipyleaflet.Map(center=location, zoom=zoom)

    # Get our starting and ending nodes
    start_node = nodes_frame.loc[route[0]]
    end_node = nodes_frame.loc[route[-1]]
    start_xy = [start_node.y, start_node.x]
    end_xy = [end_node.y,end_node.x]

    # Add both nodes as markers
    m.add_layer(ipyleaflet.Marker(location=start_xy, draggable=False))
    m.add_layer(ipyleaflet.Marker(location=end_xy, draggable=False))

    # Add edges of route
    for u,v in zip(route[0:],route[1:]):
        try:
            x,y = (ways_frame.query(f'u == {u} and v == {v}').to_dict('list')['geometry'])[0].coords.xy
        except:
            x,y = (ways_frame.query(f'u == {v} and v == {u}').to_dict('list')['geometry'])[0].coords.xy

        points = map(list, [*zip([*y],[*x])])
        ant_path = ipyleaflet.AntPath(
            locations = [*points], 
            dash_array=[1, 10],
            delay=1000,
            color='#7590ba',
            pulse_color='#3f6fba'
        )
        m.add_layer(ant_path)

    return m
draw_route(graph, shortest_route)

Retrieve buildings from named place

Just like our graph above, we can also retrieve all the building footprints of a named place.

# Retrieve the building footprint, project it to match our previous graph, and plot it.
buildings = osmnx.geometries.geometries_from_address(place_name, tags={'building':True}, dist=500)
buildings = buildings.to_crs(edges.crs)
osmnx.plot_footprints(buildings)
../_images/GettingStarted_45_0.png
(<Figure size 576x576 with 1 Axes>, <AxesSubplot:>)

Now that we have the building footprints for the campus, let’s plot that shortest route again.

First, get the nodes from the shortest route, create a geometry from it, and then visualize building footprints, street network, and shortest route all on one plot.

from shapely.geometry import LineString

# Nodes of our shortest path
route_nodes = nodes.loc[shortest_route]

# Convert the nodes into a line geometry
route_line = LineString(list(route_nodes.geometry.values))
# Create a GeoDataFrame from the line
route_geom = geopandas.GeoDataFrame([[route_line]], geometry='geometry', crs=edges.crs, columns=['geometry'])

# Plot edges and nodes
ax = edges.plot(linewidth=0.75, color='gray', figsize=(15,15))
ax = nodes.plot(ax=ax, markersize=2, color='gray')

# Add building footprints
ax = buildings.plot(ax=ax, facecolor='khaki', alpha=0.7)

# Add the shortest route
ax = route_geom.plot(ax=ax, linewidth=2, linestyle='--', color='red')

# Highlight the starting and ending nodes
ax = od_nodes.plot(ax=ax, markersize=100, color='green')
../_images/GettingStarted_48_0.png

Calculating route length

The length of the shortest route can also be determined as follows:

route_geom.length[0]
/tmp/ipykernel_11410/990676892.py:1: UserWarning: Geometry is in a geographic CRS. Results from 'length' are likely incorrect. Use 'GeoSeries.to_crs()' to re-project geometries to a projected CRS before this operation.

  route_geom.length[0]
0.011196459643642465

Notice the warning above? Our route line is projected in EPSG:4326, which uses degrees as a unit of measurement. That’s why it gives the length of our route as 0.011196459643657595 (approximately 1.23 km). This Geographic CRS isn’t accurate for geometric calculations, so we’ll need to reproject the route. Let’s project the line into EPSG:3857, which uses meters as a unit of measurement.

route_geom = route_geom.to_crs(epsg=3857)
route_geom.length[0]
1407.21083766088

Querying points of interest (POIs) using Overpass API

Overpass API is OSM’s querying API. It is incredibly powerful in that it can very quickly return queried features, and allows for selection of location, tags, proximity, and more.

Let’s query for restaurants near the University of Toronto.

import overpass

api = overpass.API()

# We're looking for restaurants within 1000m of a given point
overpass_query = """
(node["amenity"="restaurant"](around:1000,43.66, -79.39);
 way["amenity"="restaurant"](around:1000,43.66, -79.39);
 rel["amenity"="restaurant"](around:1000,43.66, -79.39);
);
out center;
"""

restaurants = api.get(overpass_query)

The example above uses the overpass package, which by default returns results in geojson format. See the overpass documentation for more information.

Next, let’s extract some data about each restaurant, and then plot all of them on a map. This time, we’ll use a plotly ScatterMapBox, which uses tiles from MapBox. You can refer to plotly’s documentation here. You can hover over each point and see the name of the restaurant (POI) that was plotted.

import plotly.graph_objects as obj

# Extract the lon, lat and name of each restaurant:
coords = []
text = []
for elem in restaurants['features']:
    latlon = elem['geometry']['coordinates']
    if latlon == []: continue
    coords.append(latlon)
    if 'name'  not in elem['properties']:
        text.append('NONAME')
    else:
        text.append(elem['properties']['name'])
# Convert into a dictionary for plotly
restaurant_dict = dict(type='scattermapbox',
                   lat=[x[1] for x in coords], 
                   lon=[x[0] for x in coords],
                   mode='markers',
                   text=text,
                   marker=dict(size=8, color='blue'),
                   hoverinfo='text',    
                   showlegend=False)


# plotting restaurants' locations around University of Toronto

center=(43.662643, -79.395689) # UofT main building

fig = obj.Figure(obj.Scattermapbox(restaurant_dict))

# defining plot layout
fig.update_layout(mapbox_style="stamen-terrain")
fig.update_layout(margin={"r":0,"t":0,"l":0,"b":0}, mapbox = {'center': {'lat': center[0], 'lon': center[1]}, 'zoom': 13})
fig.show()

Summary

Hopefully by now you’ve familiarized yourself with the various libraries and tools that are available to use in python for geospatial data manipulation and analysis. While there is plenty more to learn outside of these examples, the content covered in this section will give you a basic understanding of the possible operations that can be done on geospatial data.

For more details and library-specific best practices for each of the libraries used here, consult the respective packages’ API documentation.